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ABSTRACT 

Paneling  methods  are  approximate  techniques  for  solving  flow  pro¬ 
blems  over  wings  and  bodies.  Vortex  panels  are  used  to  model  flow  over 
wings  and  other  lifting  surfaces.  The  author  develops  a  triangular 
vortex  panel  having  a  vorticity  distribution  that  can  vary  in  magnitude 
and  direction.  This  panel  is  used  to  predict  the  pressure  distribution 
on  a  rectangular  and  a  sweptback  wing  in  subsonic  flow.  Lift  distributions 
obtained  compare  favorably  to  Anderson's  solution  and  wind  tunnel  results. 
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DEVELOPMENT  AND  APPLICATION 
OF  A  SUBSONIC  TRIANGULAR 
VORTEX  PANEL 


I .  Introduction 

Background 

Paneling  methods  are  approximate  techniques  for  solving  linearized 
subsonic  and  supersonic  potential  flow  problems  over  wings  and  bodies. 
Panels  in  use  today  incorporate  singularity  d i st ribut ions v of  the 
source,  vortex,  and  doublet  types.  Source  panels  are  often  used  to 
model  bodies  and  other  non-lifting  surfaces  and  to  model  thickness 
of  lifting  surfaces.  Vortex  panels  are  used  to  model  either  lifting  or 
non-lifting  surfaces.  In  a  typical  problem,  the  airplane  is  represented 
by  a  finite  number  of  panels.  Each  panel  is  a  singularity  distribution 
of  unknown  strength  that  models  some  part  of  the  aerodynamic  surface. 
These  unknown  strengths  are  determined  by  applying  the  flow  tangency 
boundary  conditions  at  or  near  the  aerodynamic  surfaces.  Once  the 
strengths’  are  known,  the  perturbation  velocities  can  be  computed.  These 
are  substituted  into  the  Bernoulli  equation  to  obtain  the  pressure 
distribution  and  the  corresponding  aerodynamic  forces  and  moments. 


Problem  Statement 

In  many  current  paneling  methods,  the  orientation  of  the  vorticity 
vector  is  fixed  on  the  panel  (Ref  6  and  Ref  8).  This  leads  to 
unacceptable  errors  in  some  cases.  The  purpose  of  this  study  is  to 
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derive  a  subsonic  triangular  panel  having  a  vorticity  distribution  that 
can  vary  both  in  magnitude  and  direction.  A  panel  system  is  assembled 
and  used  to  predict  the  pressure  distribution  on  a  planform.  Since  the 
computations  involved  are  too  lengthy  to  be  performed  by  hand,  a 
computer  code  was  developed  to  apply  the  panels  to  planar  wings. 

Approach 

The  subsonic  triangular  panel  is  derived  in  Section  II.  The 
Biot-Savart  Law  is  used  to  compute  the  induced  velocity  at  a  point  due 
to  an  assumed  bilinear  vorticity  distribution  on  the  panel.  This 
expression  formulates  the  induced  velocity  in  terms  of  panel  geometry 
and  unknown  corner  point  vorticity. 

Section  III  presents  methodology  for  panel  system  assembly. 
Numbering  schemes  are  developed  for  panels  and  the  unknown  corner  point 
vorticities.  Planform  boundary  conditions  are  applied  which  reduce  the 
number  of  unknowns  and  the  remaining  unknowns  are  solved  for  by 
formulating  a  linear  system  of  equations.  This  system  consists  of 
control  point  equations  (one  per  panel)  and  a  number  of  edge  continuity 
conditions.  The  linearized  form  of  the  flow  tangency  boundary  condi¬ 
tion  is  then  used  to  effect  the  solution.  Once  the  corner  vorticities 
are  known,  the  vorticity  at  any  point  on  the  planform  can  be  obtained. 
Finally,  induced  velocities  and  corresponding  pressures  can  be  calcu¬ 
lated  from  the  known  vorticity. 


Section  IV  presents  program  predictions  for  a  rectangular  wing 
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which  are  compared  against  Anderson's  solution  (Ref  1 :9~ 16) .  Pre¬ 
dictions  are  also  presented  for  a  swept  untapered  wing  which  are 
compared  against  wind  tunnel  tests  (Ref  4:92). 

Section  V  concludes  the  report  and  makes  recommendations  for  the 
improvement  of  the  aerodynamic  model. 
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I 1 .  Panel  Derivation 

This  section  presents  the  development  of  the  subsonic  triangular 
panel.  The  goal  is  to  derive  an  expression  for  the  induced  velocity  at 
an  arbitrarily  chosen  point  in  the  xy  plane  due  to  an  assumed  bilinear 
vorticity  distribution  on  the  panel. 

Geome  t  ry 

The  first  step  in  panel  development  is  the  definition  of  panel 
geometry.  Initially  assume  the  panel  is  a  trapezoid  lying  in  the  xy- 
plane  and  having  two  edges  parallel  to  the  x-axis.  It  is  then  sub¬ 
divided  into  two  triangles  having  a  common  side  that  joints  the  upper 
left  hand  corner  to  the  lower  right  hand  corner.  The  panel  is  oriented 
so  that  the  root  and  tip  chords  lie  parallel  to  the  free  stream  flow 
direction  at  CX  ■  0.  Figure  1  depicts  the  panel  geometry,  corner 
point  numbering  scheme,  and  coordinate  system. 

Singularity  Strength  Distribution 

The  general  form  of  the  singularity  strength  distribution  on  the 
panel  will  be 

Cx)(x»  y)  ■  $(x,  y)T+  y(x,  y)"j  (2*D 

where  §  and  'Y  are  continuous  functions  of  x  and  y.  For  the  purpose 
of  this  study,  §  and  are  assumed  to  have  the  "bilinear"  forms 

(2.2) 


k 


6(x.  y)  -  ?  +  Dx  +  Ey 
7U.  y)  •  A  +  Bx  +  Cy 


(2.3) 
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where  the  coefficients  A,  B,  C,  D,  F,  D  and  E  are  yet  to  be  determined 
constants.  These  constants  will  be  expressed  as  functions  of  the  panel 
geometry  and  the  unknown  singularity  strengths  at  the  corner  points. 

Consequence  of  the  Helmholtz  Condition 

The  vorticity  distribution  Cl)(x,  y)  is  required  to  satisfy  the 
Helmholtz  condition  that  vorticity  must  be  preserved  in  the  fluid.  Thus, 


v<6J  + 

0 

l 

(2.4) 

on  the  panel .  If  (§ 

1  and  have 

the  forms  (2.2)  and  (2.3)  , 

respec- 

tively,  -tfcen 

9(F  +  Dx  +  Ey)/  0x  + 

0(A  4  Bx  4  Cy)/  Qy  -  0 

(2.5) 

which  impl ies 

D-  =-  C  (2.6) 

such  that 

fi(x,  y)  -  F  -  Cx  4  Ey  (2>7) 

The  formulation  could  be  continued  in  terms  of  the  bilinear 
coefficients.  However,  it  is  conventional  to  express  them  in  terms  of 
panel  corner  point  vorticities. 

Bilinear  Coefficients  in  Terms  of  Corner  Vorticity 

In  the  ensuing  discussion,  the  subscript  L  refers  to  the  leading 
triangle  and  the  subscript  T  refers  to  the  trailing  triangle. 

The  'Y  component  of  vorticity  (eq  (2.3))  is  assigned  the  unknown 
values  %  %.  %  at  the  corners  (x^  ,  yj),  (xj,  y2) ,  and  (xj,,  y2) 
of  the  leading  triangle.  This  leads  to  the  following  system  of  three 


6 


ASO-TR-81-5005 


equations  in  AL,  B|_  and  CL : 

XmAL  +  h*i  +  ci7i  (2.8) 

^3  "  \  +  *LX3  +  °Ly2  <2*9> 

X  “  *L  *  *tx4  +  °Ly2  <2*1°) 

The  system  has  the  solution: 

\  -  72  \^72  “  yl>  +  1<ylx4  “  xly2)/(y2  "  yl)5 
^/(*3  -  x4)  +  [  (Xjjj  -  yix3^/(y2  “  yl)J 

ty(*3  -  «4)  (2.1D 

Bl  -  (  73  -  X>/<*3  -  x4)  (2.12) 

-  ^/(y2  "  yl5  +  l(xl  “  x4)/(y2  “  yl)5  %/ 

(*3  -  x4)  +  l(*3  -  Xj)/(y2  -  yx)l  %/ 

(x3  -  x4)  (2.13) 

Since  DL  =  ”  C|_>  only  two  corner  conditions  may  be  used  to  solve  for 
and  El  in  the  §  component  eq  (2.7).  Assigning  the  values  §,  and  63 
at  the  corners  (xj  ,  y^)  and  (x^,  y2)  leads  to 

4  +  °LX1  "  Vl  +  FL  (2-14) 

S3  +  °LX3  *  \y2  +  FL  (2*15> 
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*l  "  (  $1  -  53)/(71  -  y2)  +  [(Xj  -  *3)/ 

(7X  “  ypiCj.  (2 .16) 

FI  "  (yl  63  "  ^2  5^/(7!  -  72)  +  [(7^3  - 

Vp/Cyj.  *  y2)1CL  (2.17) 

In  a  similar  way,  the  trailing  triangle  coefficient  equations  are 

obtained  as  follows: 

“  yl  Ti^i  ■  72)  +  U?2*l  ~  X4yl)/(yl  ■  J2)] 

72/(x2  -  *x)  +  I  (x47x  -  72x2)]  Oy 
(x2  -  Xj)  (2.18) 

Bt  -  (  72  -  7i>/(*2  -  x2)  (2.19) 

S  *  -  ^4/(73  -  72>  +  I(x4  -  x1)/(y1  -  y2)]  y2/ 

(x2  -  Xj)  +  I(x2  -  x4)/(y2  -  y2)]  7^/ 

(x2  -  *j)  (2.20) 

*T  “  (  62  -  ^5/(7!  -  72>  +  l(x2  -  x4)/ 

(7i  -  72)Kt  (2.21) 

By  ■  (72  5  4  -  72  5  2)Z(y2  ~  y2^  +  ^ylx4  ** 

x272)/(71  "  72>3Ct  (2‘22) 
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where  Ey  and  Fy  have  been  expressed  in  terms  of  (as  functions  of)  (jJ 
components  at  corners  2  and  k. 

Mathematical  Continuity 

The  bilinear  vorticity  distribution  is  continuous  on  each  triangular 
region.  In  addition,  the  component  has  been  made  to  be  continuous 
throughout  the  planform  by  the  representation  in  terms  of  corner  values. 
This  can  be  demonstrated  as  follows.  Let  J"*  be  the  boundary  shared  by 
any  two  adjacent  triangles.  Then  J"1  is  a  straight  line  segment  and  can 
be  described  by  a  linear  expression  (i.e.,  y  in  terms  of  x  or  x  in  terms 
of  y).  The 'y  distribution  on  each  of  the  adjacent  triangles  will 
degenerate  to  a  linear  function  of  a  single  variable  upon  substitution 
of  this  expression.  Both  functions  assume  the  same  /Y  values  at  the 
endpoints  of  £"*  .  Since  only  two  points  are  needed  to  determine  a 
straight  line  or  a  linear  form,  we  have  ~Y  matching  identically  on 

r . 

The  (5  component  has  breaks  in  continuity  throughout  the  planform. 
This  is  a  consequence  of  applying  the  Helmholtz  condition  (eq  (2.1*)) 
which  eliminated  the  constant  D  and  expressing  the  remaining  two 
unknowns  in  terms  of  two  ($  corner  values,  out  of  a  possible  three. 

The  (5  component  is  continuous  on  panel  leading  and  trailing  edges 
since  it  is  on  these  edges  that  common  (5  values  are  assumed  at  the 
endpoints.  The  (5  component  is  discontinuous  on  panel  diagonals. 

Application  of  the  Biot-Savart  Law 

Let  (jJ(x)  be  a  vorticity  distribution  defined  on  a  finite  region 
R  in  the  xy  plane.  Let  $  be  a  fixed  point  (control  point)  in  the  plane. 
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The  velocity  induced  at  s  due  to  the  distribution  on  R  is  given  by  the 
Biot-Savart  Law  (Ref  *5 : S26-l>28^ : 


•//' 


47T?<*)  -  l  G(t)  X  (S  -  X)]/|t  - 


(2.23) 


Suppose  the  control  point  s  is  located  at  the  origin  of  the 
coordinate  system.  Note  that  this  can  be  done  by  performing  a  simple 
translation  of  the  plane.  Then, 

•  -  (o,  0,  0)  (2.24) 


•  -  S  -  (0,  0.  0)  -  (x,  y,  0)  -  (-x,  -J,  0) 

Assuming  the  distribution  has  the  form  (2.1), 

Qx  cr  -  *)  -  )k 


(2.25) 


(2.26) 


where  k  is  the  unit  vector  normal  to  the  xy-plane.  Substituting  the 
expressions  (2.25)  and  (2.26)  into  eq  (2.23)  yields: 


47Tw  (o,  0,  0)  mJJ (xY  -  J  ft  )/(x2  +  y2)3/2dR  (2.27) 


R 


Where  w  is  the  normal  velocity  component  induced  at  the  origin. 
Substituting  the  expressions  for  ^(2.3)  and  (5(2.7)  into  eq  (2.27) 
yields : 


n'-fl 


(Ax  +  Bx  +  2Cxy  -  Ey2  -  Fy)/ 


(X2  +  T2)V2dR 


(2.26) 


where  R  is  taken  as  the  region  defined  by  a  trapezoidal  panel. 


Let  R^  and  R-j.  be  the  subregions  of  R  which  correspond  to  the  lead 
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ing  and  trailing  triangles.  The  coefficients  A,  B,  C,  F,  and  E  remain 
constant  on  each  subregion  and  eq  (2.28)  is  rewritten  as: 

*71-..  vi  +  v*  +  v\  *  vl  +  IC^  ♦ 

2Vs  -  V*  -  vl  -  Vs  -  Vs  «-2» 

where 

X3*  T  -  //  +  72)V2W 

V  T 

12>  T  “  //  (x2/(*2  ♦  y2)3/2>dR 

\  T 

Ij’  T  -  JJ  (*y/(*2  +  y2)3/2)dR 

*L,  T 

Ii*  T  -  JJ  <y/(*2  +  72)3/2)dR 

*L,  T 

25*  T  -  //  fr2'<*2  *  y2>3/2>« 

V  T 

Evaluation  of  these  integrals  may  be  found  in  Appendix  A.  Substitution 
of  the  expressions  for  A^  through  Ey  and  collecting  coefficients  of  the 
unknowns  ( 6,,  §2, 63, 64.  7,,  y2,  7 and  O4)  yields,  after 
considerable  algebraic  manipulation, 


(2.30) 


(2.31) 


(2.32) 


(2.33) 


(2.36) 
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4  7T  *  -  t(y  Jl4-  ify(yx  -  y2)]  6X  + 

C(y2iJ  -  i*)/^  -  y2)l  82  + 

«I^  -  y1lj)/(y1  -  y2)]  63  + 

UI5  -  y^/c*!  -  y2>]  S4  + 

I^Vl  ”  y2x2^Xl  ~  ^x2  “  xl?y2Xl 

(y2  -  y2)i2  (*2  “  + 

(x2  -  »4)kt>/((x2  -  »2) (yx  -  y2))l  'Y1  + 

I«72*l  -  *47x)lJ  ♦  <JX  -  7j»2  + 

(*4  -  x^J/Kyj  -  y2)(*2  -  *!»]  y2  + 

I((yl*4  '  Xly2)Xl  +  (72  -  yl)X2  + 

<*3  -  *A)\)/((y2  “  Ti> -  V51  %  + 

*(*xly2  "  yl*3^Xl  ”  ^x3  "  x4^ylXl  " 

(y2  “  y2)lj  +  (*3  -  Xj)\  +  <*3  -  X4>V7 

((y2  -  yx)(*3  -  *4>)]  ^4 


(2.35) 
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where 


*L  "  ^3  ~  ^X1  “  x3^^1  ”  y2^*5  “ 

[(71*3  -  *iy2>/(yi  -  72*11A  (2.36) 

*1  -  2I3  -  “*2  -  V/O-i  -  72)>iJ  - 

1(7^4  -  *2y2)/(yl  "  y2^X4  (2.37) 


Summary 

Expression  (2.35)  is  the  normal  velocity  induced  at  the  origin  of 
the  xy  plane  by  a  trapezoidal  vorticity  panel.  This  velocity  is  due  to  a 
bilinear  vorticity  distribution  which  satisfies  the  Helmholtz  condition 
eq  (2.U). 
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III.  Panel  Assembl' 


This  section  presents  the  panel  assembly  procedures  needed  to 
model  a  wing.  The  goal  is  to  develop  the  methodology  required  to 
predict  the  pressure  distribution  and  associated  forces  and  moments  on 
a  wing. 


Panel  Numbering 

Figure  2  illustrates  a  paneling  arrangement  and  associated  number¬ 
ing  scheme  for  a  16  panel  wing.  The  panels  are  numbered  consecutively 
in  the  chordwise  direction  starting  with  the  inboard  leading  edge 
panels  and  terminating  with  the  outboard  trailing  edge  panels. 


Number  of  Unknowns,  Boundary  Conditions  and  Numberin 


Let  M  be  the  number  of  chordwise  panels  and  N  be  the  number  of 
spanwise  panels.  These  are  defined  using  M  +  I  chordwise  cuts  and 
N  +  1  spanwise  cuts.  Each  intersection  determines  a  panel  corner 
point.  Since  there  are  two  unknown  components  at  each  corner  point, 
the  total  number  of  unknowns  is  given  by: 

2(M  +  1)  (N  +  1)  (3-D 


Boundary  Conditions 

Two  boundary  conditions  are  imposed  on  the  wing  panel  system. 
These  reduce  the  number  of  unknowns  and  improve  the  physical  modeling 
of  the  flow  field. 

The  Kutta  condition  (Ref  5:390-399) 
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Figure  2.  A  Fueling  Arrugenent  end  Associated 
Numbering  Scheme  for  a  16  Panel  Ving 
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T(x,  y)  =  0 


(3-2) 


is  imposed  at  all  corner  points  on  the  wing  trailing  and  tip  edges. 
Since  the  'Y  component  of  vorticity  eq  (2.2)  is  both  linear  and 
continuous  on  these  edges,  the  Kutta  condition  is  satisfied 
identically.  The  Kutta  condition  reduced  the  total  number  of 
unknowns  by: 

M  +  N  +  1  (3.3) 

(NOTE:  The  point  defining  the  intersection  of  the  trailing  and  tip 
edges  is  common  to  both.) 

A  historically  acceptable  boundary  condition  (an  outgrowth  of 
Prandtl's  lifting-line  theory  (Ref  5=535“567))  is  for  the  vorticity 
vector  to  lie  tangent  to  the  wing  leading  edge.  This  boundary  condi¬ 
tion  initially  orients  the  vorticity  vector  so  that  a  positive  circu¬ 
lation  is  produced.  The  boundary  condition  is  imposed  at  all  leading 
edge  corner  points  and  can  be  written: 


or 


7/8  -  A 


5  *  7/ A 


(3-M 


(3.5) 


where  A.  is  the  leading  edge  slope  at  the  corner  point.  It  reduces 
the  total  number  of  unknowns  by: 


N  +  I 

and  the  total  number  of  unknowns  becomes: 


(3-6) 
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2MN  +  M  (3.7) 

after  imposing  the  two  boundary  condition  equations  (3.2)  and  (3*5). 

Since  the  wing  is  symmetric  about  the  x-axis,  it  may  seem  logical 
to  impose  a  boundary  condition  at  the  wing  root  chord.  However,  setting 

5=0  (3.8) 

at  the  centerline  is  redundant  for  rectangular  wings  and  leads  to  an 
ill-conditioned  system  once  planform  symmetry  is  considered. 

Unknown  Nuirfcerinq  Scheme 

Figure  3  illustrates  the  unknown  numbering  scheme  for  a  9  panel 
wing  with  applied  boundary  conditions.  The  paneling  arrangement  of 
Figure  3  is  chosen  because  it  represents  the  smallest  number  of  panels 
needed  to  illustrate  interior  panels  and  panels  having  boundary 
conditions.  Let  i  be  the  panel  number  of  an  interior  panel.  Then 
the  following  numbers  (in  terms  of  i,  M,  and  N)  are  assigned  to  the 
eight  unknown  corner  vorticity  components: 


< 6 »!> 

1-1 

(3.9) 

<S2> 

1 

(3.10) 

(63) 

i  +  M  -  1 

(3.11) 

(64) 

i  +  M 

(3.12) 
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The  unknown  corner  5  s  are  denoted 
by  and  the  7a  by  "■ 


X 

Figure  3.  Unknown  Numbering  Scheme  for  a  9  Panel  Wing 
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<7i> 

MN  4 

M  4  1 

(3.13) 

<7Z) 

MN  4 

M  4  1  4  1 

(3.14) 

<r3) 

MR  4 

2M  4  1 

(3.15) 

<X> 

MN  4 

2M  4  i  4  X 

(3.16) 

■  the 

Corner 

Vortici ties 

The  total  number  of  unknowns  (after  the  boundary  conditions  are 
applied)  is  given  by  eq  (3-7)  which  also  specifies  the  number  of  condi¬ 
tions  needed  to  solve  for  the  corner  vorticities.  Two  types  of  condi¬ 
tions  will  be  used;  control  point  conditions  and  edge  or  point 
continuity  conditions. 


Control  Point  Equations 

Control  point  equations  are  obtained  using  eq  (2-35).  The 
velocity  component  is  computed  for  one  point  (control  point)  on  each 
panel  comprising  the  wing. 


Let  (x j ,  yj)  be  the  control  point  on  panel  i.  To  obtain  the 
contribution  to  Wj  due  to  panel  j,  express  the  coordinates  of  panel  j 
in  a  coordinate  system  with  (xj ,  yj)  at  the  origin.  This  is  done  by 
performing  a  translation  in  the  z  ■  o  plane: 


x*  ■  x  -  x. 


Equation 


(2.35) 


(3.17) 

y '  -  y  -  y  j 

is  then  applied  with  appropriate  boundary  conditions. 
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Planform  symmetry  is  included  by  reflecting  either  the  panel  or  con¬ 
trol  point  about  the  x-axis  and  applying  eq  (2.35)-  Reflecting  the 
control  point  is  less  complicated  from  a  programming  viewpoint. 

The  above  process  is  repeated  for  each  pane)  on  the  wing.  After 
all  the  contributions  to  w.  have  been  calculated,  it  can  be  written  as: 

2MN  +  M 

w±  ’  Z  (3,l8> 

J  -  1 

where  0 •  is  the  column  vector  of  unknown  corner  vorticities  and  the 
coefficients  Ajj  are  functions  of  panel  geometry.  The  are  numbered 
using  the  system  given  by  eqs  (3.9)  through  (3-16).  One  control  point 
equation  (3*18)  is  obtained  for  each  panel  on  the  wing  and  together 
they  comprise  MN  conditions. 

Edge  Continuity  Conditions 

The  (5  component  of  vorticity  eq  (2.7)  is  discontinuous  across 
the  panel  diagonal  (See  discussion  in  Section  II.).  This  can  be 
partially  remedied  by  specifying  a  point  continuity  condition  at  the 
panel  lower  right  hand  corner.  This  condition  is: 

§l<v  y 2)  “  84  (3.»> 

which  becomes 


&  "63  +  I(x3  ~  *4>/<y2  “  yl)]  \  " 

f(xx  -  *4>/(y2  -  7l)J  y3  *  IC*3  -  *!>/ 

<72  “  *!>]  H  "  0  (3.20) 
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after  substituting  of  (x^,  yj)  and  the  expressions  for  CL  (2.13),  EL 
(2.16),  and  FL  (2.17)  into  eq  (2.7).  One  cond-tion  eq  (3*20)  is 
formulated  for  each  panel  on  the  wing  for  a  total  of  MN  condition. 

Note  that  §  is  still  not  continuous  across  the  panel  diagonal  due 
to  the  remaining  discontinuity  at  the  upper  left  hand  corner.  Also 
§is  not  continuous  across  panel  side  edges. 

The  edge  continuity  conditions  and  control  point  equations  total 
2MN  conditions.  M  additional  conditions  can  be  obtained  by  specifying 
an  edge  continuity  condition  for  at  the  upper  left  hand  corner  of 
each  panel  along  the  centerline.  This  condition  is 


yi?  "  (3-2i ) 

which  becomes 

62  -  s1  -  (<« 2  -  *j)/(Tj  -  y2))  r4  + 

-  BjJ/fyj  -  72>]  r2  +  t(*2  -  *4)/ 

(yl  "  y2^Ti  *  0  (3.22) 

after  substitution  of  (#j ,  yj)  and  the  expressions  for  Cj  (2.20),  Ej 
(2.21),  and  Fj  (2.22)  into  eq  (2.7)-  This  choice  is  based  on  trial  and 
error,  the  additional  (5  continuity  on  the  centerline  having  the  effect 
of  minimizing  vorticity  oscillations. 

Compress ibi 1 i ty 

Compressibility  effects  are  accounted  for  by  using  the  Prandtl- 
Glauert  transformation  (Ref  2: 124-127) : 
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*  “  -  K2 


(3.23) 


The  transformation  is  applied  to  all  x  coordinates  which  are  used  in 


either  the  control  point  equations  or  edge  continuity  conditions. 


Matrix  Formulation 


The  control  point  equations  and  edge  continuity  conditions  are 
2  MN  +  M  equations  in  the  unknowns,  Qj.  This  system-has  the  matrix 


formulation : 


2MN  +  M 


(ZMN  +  M)x(2MN  +  M) 


(3-24) 


The  first  MN  rows  of  A.,  are  the  coefficients  for  the  control  point 

'J 


equations  and  are  all  nonzero.  The  last  MN  +  M  rows  of  A.^  correspond 


to  the  homogeneous  edge  continuity  conditions  and  have  no  more  than 


five  nonzero  entries  per  row. 


Sol ut i on 


The  linearized  form  of  the  flow  tangency  boundary  condition  is 


(Ref  5:495) 


w/vc5 dc/d*  -a 


(3.25) 


Where  ais  the  wing  angle  of  attack,  dc/dx  is  the  local  camber  slope 


and  !/„  is  the  free  stream  velocity.  This  expression  is  substituted 


for  each  w.  in  eq  3-24  where  dc./dx  is  the  panel  slope  at  control 


point  i.  In  matrix  notation, 
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Mb) 


v^d^/dx  -or) 

O 


(3.26) 


fyvco  -  [*u] 


-1  dct/dx -Of 

o  i 


(3.27) 


Forces  and  Moments 

Once  the  ft/V  are  obtained  for  a  given  Q(  and  camber  slope 

distribution,  eqs  (2.3)  and  (2.7)  can  be  used  to  calculate  the 
vorticity  strength  at  any  point  on  the  planform.  The  surface  perturba¬ 
tion  velocities  in  terms  of  local  vorticity  are  (Ref  5=508): 


o/VCO  - 

±  T/2 

(3-28) 

T/Too- 

+  6/2 

(3-29) 

The  upper  sign  corresponds  to  the  upper  wing  surface  and  visa-versa. 

Pressure  coefficients  are  obtained  from  the  perturbation  velocities 
by  either  using  the  exact  isentropic  expression  (Ref  3:167): 


Cp  -  2[(1  +  (7-  1)M^/2(1  -  ((V^-F  u)2  ♦ 

▼2  ♦  v2)  y^y- x)  -  u/^^) 

or  the  second  order  approximation  (Ref  3=167): 

Cp  •  *  t2u/V00+  (1  -  >4))u2/V00+  (V2  +  v2)/V2  ] 


(3-30) 


(3.31) 
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which  is  adequate  for  two-dimensional  and  planar  flows.  These  coeffi¬ 
cients  are  integrated  along  chord  lines  to  obtain  local  lift  and  moment 

coefficients.  The  appropriate  expressions  are 

’*TE 

c,  ■  (1/c)  I  (C_  -  C_  )d* 


r  xe 

o 


0.32) 


c«  -  <1/e' 


>[  <VVS 

\e 


dx 


(3-33) 


where  the  subscripts  1  and  u  refer  to  the  lower  and  upper  wing 
surfaces . 
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IV.  Results 

A  FORTRAN  Code  "WING"  was  developed  to  analyze  planform  flow  using 
the  methodology  discussed  in  Sections  II  and  III.  Most  of  the  program¬ 
ming  techniques  used  in  WING  were  previously  developed  by  the  author 
and  can  be  found  in  Reference  8. 

Program  WING  was  exercised  for  a  variety  of  four  point  wings, 
taper  ratios,  and  sweep  angles.  This  section  presents  results  for  two 
of  these  wi ngs . 

Two  general  observations  are  made  first.  One,  control  point 
location  is  the  major  factor  controlling  bounded  numerical  oscillations 
of  the  vorticity  vector  as  it  is  in  many  current  paneling  routines 
(ex.  Refs  6  and  8).  Oscillations  are  very  common  if  the  control  point 
is  located  anywhere  on  the  leading  triangle.  Fewer  oscillations  occur 
if  the  control  point  is  located  on  the  trailing  triangle  with  . 1  *  CY  * 
.5  and  .A  *  CX  «  .9.  Secondly,  the  program  shows  the  best  results  when 
uniform  spanwise  paneling  is  used.  Non-uniform  spanwise  paneling  tends 
to  cause  oscillations  in  the  vorticity  vector.  However,  non-uniform 
chordwise  paneling  seems  to  have  little  effect  on  solution  stability. 
The  best  total  match  (with  other  known  solutions)  occurs  at 
approximately  CY  *  .15  and  CX  ■  .75. 

Rectangular  Wing 

The  first  case  examined  is  a  rectangular  wing;  AR  ■  8, Of  *  5° 
and  Mqq  *  .1.  The  wing  is  modeled  using  11  uniformly  spaced  span 
stations  and  6  non-uni formly  spaced  chord  stations  (0.,  .1,  .3,  .6, 


25 


4f« 


V- 


, 


I 


ASD-TR-81 "5005 

.8,  and  1.)  which  define  50  panels.  Figure  4  shows  the  Cl  distribution 
predicted  by  WING  and  Anderson  (Ref  1:9-16).  Figure  5  shows  the  span- 
wise  distribution  of  center  of  pressure.  The  Xcp  stays  close  to  the 
classical  25%  chord  location  except  in  the  outboard  region  where  it 
travels  forward  as  expected. 

Figure  6  is  the  y strength  at  the  root  chord.  The  calculated 
solution  compares  favorably  with  the  exact  2-D  flat  plate  solution 
(Ref  5:515): 

T(x)  =  2d(c  -  x)/(cx  -  x2)  (5.1) 

where  c  is  the  chord  length  andQf  is  measured  in  radians. 

Swept  Wing 

The  second  case  examined  is  a  swept  untapered  wing;  AR  =  4.5, 

A*  40 °,Ct=  5°,  and  Mqq  =  The  wing  modeling  is  the  same  as  the 


rectangular  wing.  Figure  7  shows  the  Cl  distribution  predicted  by 
WING  and  by  wind  tunnel  tests  (Ref  4:92).  Figure  8  shows  the  x  center 
of  pressure  versus  span  station.  Again,  good  agreement  is  obtained  for 
both  the  C.  distribution  and  x  travel. 

L  r.o 
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V.  Conclusions  and  Recommendations 


Concl us  ions 

A  paneling  technique  which  allows  the  vorticity  vector  to  change 
direction  has  been  demonstrated.  The  proper  choice  of  control  point 
(see  section  IV)  will  guarantee  a  solution  free  of  numerical  oscil- 
lations  in  the  vorticity  vector.  This  method  has  been  implemented 
on  a  computer  and  introduces  no  new  complexities  to  an  experienced 
programmer.  Existing  mesh  generators  and  other  aerodynamic  modules 
were  incorporated  into  this  technique  (Ref  8).  Computer  run  times 
are  of  the  same  order  as  programs  incorporating  "fixed  direction" 
vorticity  panels  and  no  problems  involving  extensive  "run  times"  were 
encountered. 

The  method  gives  satisfactory  aerodynamic  results  for  the  cases 
examined.  However,  as  in  other  paneling  methods,  the  results  are 
control  point  sensitive  and  some  "curve  fitting"  is  involved  to  find 
the  control  point  location  that  best  matches  existing  data  or  solutions. 
Also,  the  best  "fits"  don't  shown  any  appreciable  improvement  over 
existing  paneling  methods  such  as  the  constant  strength  vorticity 
panel  program  described  in  reference  8. 


Recommendat i ons 

The  following  ideas  are  suggested  for  further  study  or  improvement 
of  the  method: 

a.  The  method  should  be  tested  for  cranked  wings. 

b.  The  panel  should  be  oriented  in  space  to  allow  method  assess¬ 
ment  for  three  dimensional  flow  situations. 
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APPENDIX  A 

Evaluation  of  the  Panel  Integrals 

The  methodology  for  evaluating  the  integrals  l*T’^  of  section  2  is 
presented  here. 

Consider  the  semi - i nf i n i te  region  shown  in  Figure  9-  This  strip 
is  bound  by  the  lines  y  *  yQ,  y  *  yj,  and  the  line  segment  connecting 
(xq,  yQ)  to  (x, ,  yj).  The  equation  for  the  line  segment  is: 


y  •  Mr  +  b 

(A.  1 ) 

where 

M° 

1 

& 

i 

H 

£ 

■ 

X 

(a. 2) 

and 

*>  *  7.  -  *  M 

0  0 

(A.3) 

Define  the 

following  improper  integrals  on  the  semi 

-infinite 

rTi  ft* 

I  I  7)dxiy  (i  -  1,  - 5) 

L*OoJjo  J  (y  -  b)/* 

(a.i») 

where 

y)  ■  */<x2  +  y2j3/2 

(A.  5) 

*2<*'  7)  -  *2/(*2  +  y2)3/2 

(A.6) 
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*3<*»  y)  -  xy/(x2  +  y2)3/2 

(A.  7) 

y)  -  y/(x2  +  y2)3/2 

(A.8) 

*5<*.  y)  -  y2/(*2  +  y2)3/2 


(A.  9) 


To  illustrate  the  methodology  used  In  evaluating  the  integrals  (A. 4) , 
consider 


11a 
Ir*  CO 


L 

(y 


ty2/(x2  +  y2)3/2]dxdy 
-  b)/M 


(A.  10) 


By  Pierce's  intetral  tables  (Ref  7) , 


-5  -  Ua  f  1  [y2x/(y2V?+  yS 


L*CO  y 


(a. id 


(y  -  b)/M 


llm 


]>COy  y 


J  L  dy/VL2  +  y2  -  J  1  «y  -  b)/M)dy/ 


V«7-  b)/M)2  +  y2 


(A. 12) 


Substi tut i ng 


*  -  (y  -  b)/M 


(A. 13) 


and 


dy  ■  Kdz 


(A. 14) 
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into  the  rightmost  integral  of  (A. 12)  yields: 


?5  -  liffl 


r 


(Ldy/V^2  +  y2)-H 


1^00  y 


f 


zdz/ 


^7 


(Mz  +  b)‘ 


(A.15) 


The  integrals  (A.15)  can  be  evaluated  using  Ref.  7.  One  obtains,  after 
some  alg'ebraic  manipulation, 

*5  -  11a  L  +  ^|Ti2)/(yo  +  *^T"l2)]  + 

L+00 

twVa  +  +  ((xx  +  Hyp/ 

V1  +  +  yo  4  ((xo  4  ^0)/V1  +  m2))J  + 

tH/u  +  ^KA^TT;-  vCf77f)  (A.,W 

recal 1 ing 

yo"Msto  +  b  (A. 17) 

and 


yl-Hx1  +  b  (A.  18) 

Notice  the  limit 

Urn  L  laKyj  +  'JJ[+l2)/(70  +  *^2  +  L2))  (a.  19) 

L«*00 
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is  not  a  function  of  xQ  or  xQ  which  leads  to  the  following  observation. 
If  the  integral  is  evaluated  on  any  other  semi-inf ini te  strip  bounded 
by  the  lines  y  ■*  yQ  and  y  *  y j ,  the  limit  (A. 19)  is  invariant.  Evalua¬ 
ting  the  integral  F using  any  other  line  segment  connecting  y  *  yQ  to 
y  **  y j  and  formulating  the  difference  between  this  result  and  (A.  16) 
leads  to  cancellation  of  the  limit  (A. 19). 

The  other  four  integrals,  F|  through  F^,  can  be  evaluated  by  a 
similiar  procedure  using  Ref  7«  Each  integral  has  a  limit  term  given 


CFX)  11*  ln[(yo  +Vyf  +  l2)/(y1  +  +  L2)] 


L"*CO 


(P2)  H»Iy1l*(Vl*2  +  l\  ♦  L>  -  ypln(VL2  +  72  + 

h+OO 


(A. 20) 


(ty  u*r  jit  +  jl  -  Vl2  +  jJ) 


(P4)  ll*ln[(L  +  Vl2  +  72)/a  +  aA2  ♦  yJ)J 


L-*00 


C A .21) 


(A. 22) 


(A. 23) 


These  terms  cancel  upon  the  formulation  of  Integral  differences. 
Retaining  the  finite  terms  from  each  of  the  Integral  evaluations,  we 
can  define  "Strip  Functions"  Sj  for  the  functions  fj  and  the  points 
(x0,  yQ)  and  (xj,  y,)  by 


®it(v  v»  *i>j  •  o*/VTr?)c 


(A.  24) 
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»2I(V  *o>»  (V  V1  "  V*[V*o  +  7o  +  Xo1  “ 
Tj^lnl  ♦  yj  +  Xj]  +  (M/(l  4  H2))!!  4 


[b/(l  4  M2)3/2]6 


(A. 25) 


83I<xo,  7e),  (x1,  yx)]  -  (^/(l  4  M2))H  4 

[bM/(l  4  M2)3/2JG  . 

(A. 26) 

V(V  yo)f  (*i»  V1  •  +  V*i  +  y[)/ 

<*o  W**  +  7o)J  “  C/V^+  M2  (A. 27) 


*5I<Xo»  V*  <*1»  V1  "  &?/<*  4M2)3/2]G  - 
(M/(l  4  M2)^ 


(A. 28) 


where 

C  -  ln[ <VX1  +  yj  +  ((*!  +  +  *?))/ 

<  "V^o  +  yo  +  «*„  + 


(A. 29) 


and 


H 


(A. 30) 


Examing  Figure  1,  It  Is  obvious  that  each  of  the  panel  integrals 
(2.30)  through  (2.34)  can  be  obtained  by  formulating  the  difference  of 
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two  corresponding  strip  functions.  This  leads  to  the  following 


fundamental  results: 

T^)  »  (*3»  Tj)l  "  y^)  > 

(*4. 

l i  7j)»  (x^»  7j)]  "  Cxj i  Jj) • 

(*4.  y2)J 

where  i  ranges  from  1  to  5* 


(A.31) 


(A. 32) 
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APPENDIX  B 
Simple  Wake  Model 


The  Kutta  condition  (Ref  5:390-399)  implies 

G 0(x,  y)  =  6(x,  y)  T  (B.l) 

for  all  points  on  the  wing  trailing  edge.  Let  vortex  filaments 
emanate  from  corners  2  and  4  of  a  trailing  edge  panel  and  extend  to 
infinity  (see  Figure  10).  By  Helmholtz's  theorems,  these  filaments 
will  retain  the  strengths  §  2  an<*  64  indefinitely.  Now  assume  that 
the  strength  of  the  wake  region  between  the  two  filaments  is  given  by: 

6(*.  y)  =  (My  +  b)  T  (B.2) 

where 

M  =  ^  <5l*  62)/(y2  ■  y^  (B.3) 

b  =  d2~  y  1 6^/(y2  _  yj)l  (b*1*) 

Applying  the  Biot-Savart  Law  (Section  II),  the  effect  of  the  wake 
region  upon  the  induced  velocity  at  a  control  point  located  at  the 
origin  of  the  coordinate  system  is 

i»7Tw  *  -  M  S5[(x2,  yj),  (x^,  y2)I 

-  b  Sjj[(x2,  y,),  (x^,  y2)]  (B.5) 

The  limit  terms  (A. 23)  and  (A. 19)  for  the  "Strip  Functions"  S4  and  S5 
must  now  be  evaluated  and 

(S4)  lim  In  [(L  +  V«-2  +  yf)/(L  +^l2  +  y\)J  -  0  (B.6) 

L-*00 

($5)  Urn  L  In  [(y2  +  +  L)/(yj  +  Vy?  +  L)]  "  Y2  “  Yl  (B. 7) 

L-»00 
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Collecting  coefficients  of^  and  ,  we  obtain 
4fTw  -  [($5  -  y2Slt)/(y2  -  y1 )]  (5 2 
+  [(Y|S4  -  S5)/(y2  -  yj)]  S 4 


(B.8) 


which  is  the  induced  velocity  at  a  point  due  to  a  wake  region  emanating 
behind  an  arbitrary  trailing  edge  panel- 
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